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ABSTRACT 

In addition to fitting the data of 233 extra-solar planets with power laws, we 
construct a correlated mass-period distribution function of extrasolar planets, as 
the first time in this field. The algorithm to generate a pair of positively correlated 
beta-distributed random variables is introduced and used for the construction of 
correlated distribution functions. We investigate the mass-period correlations 
of extrasolar planets both in the linear and logarithm spaces, determine the 
confidence intervals of the correlation coefficients, and confirm that there is a 
positive mass-period correlation for the extrasolar planets. In addition to the 
paucity of massive close-in planets, which makes the main contribution on this 
correlation, there are other fine structures for the data in the mass-period plane. 

Key words: extrasolar planets, distribution functions, beta distributions, bootstrap 
method, correlation coefficient 



1. Introduction 



The observational effort has led to the discovery of more than 200 extrasolar planets 
(exoplanets), in the mass range 0.03 to 20 Jupiter Mass (Mj), with orbital periods from a 
few days to about 4000 days. Many interesting problems about the formation and evolution 
of planetary systems have therefore been studied with the information provided by these 
detected systems (Jiang & Ip 2001, Laughlin & Chambers 2001, Kinoshita & Nakai 2001, 
Gozdziewski & Maciejewski 2001, Ji et al. 2002, Ji et al. 2003, Jiang & Yeh 2004) 

Moreover, due to the growing number of detected extrasolar planets, several groups 
have been working on the statistical distributions and possible correlations. Assuming that 
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the mass and period distributions are two independent power-law functions, Tabachnik and 
Tremaine (2002) used the method of maximum likelihood to determine the best power-index. 
Though they found that the uncertainties in the mass and period distributions are coupled, 
the study of the possible mass-period correlations is beyond their scope due to the principal 
assumption of two independent power-law functions. 

In addition, Zucker and Mazeh (2002) calculated the correlation coefficient between 
mass and period for the detected data in the InP — InM space. They used Monte Carlo 
simulation to determine the p- value for testing whether the correlation is significant or not. 
They concluded that the mass-period correlation is significant. However, at the time of their 
study, the number of detected exoplanets was limited, so only 66 planets were used. 

Over the years, much more exoplanets with different properties are discovered. For 
example, more hot Jupiters are found due to the effort of transit surveys, a few newly 
detected exoplanets are moving on extremely eccentric orbits, and the exoplanets with mass 
in the order of Earth Mass are also discovered. These results have in fact brought this field 
into a completely new era. Jiang et al. (2006) did cluster analysis on 143 samples and found 
that the data grouping could be related with the dynamical processes of planetary systems. 
This approach was agreed by Marchi (2007), in which an extrasolar planet taxonomy was 
presented. 

Therefore, it is about the time to construct a new distribution function. As in Tabachnik 
and Tremaine (2002), assuming the mass and period are two independent power-law distribu- 
tions, we use 233 samples of exoplanets from exoplanet website ( ]http: / /exoplanet.eu/catalog-all.ph"p| ) 
on 6th July 2007 to construct updated distribution functions. Our samples do not include 
OGLE235-MOA53 b, 2M1207b, GQLup b, HD 187123 c, ABPic b, SCR 1845 b, SWEEPS- 
04, due to the missing of either their periods or mass. The three outliers, PSR 1257+12 b, 
HD 154345 b, PSR B1620-26 b with either extremely small mass or huge periods are also 
excluded. (Please note that the data of mass means the value of projected mass in this 
paper.) 

Moreover, because the mass and period are likely to be correlated according to the 
results of Zucker and Mazeh (2002), a new distribution function without the assumption 
that the mass and period are independent would be more satisfactory. That is, we hope to 
construct a new distribution function, in which the mass and period can be coupled. This 
was not possible until an algorithm for generating two positively correlated beta-distributed 
random variables was provided in Magnussen (2004). We therefore have to employ the beta 
distribution for this part of calculations. 

After that, in order to have a careful investigation on the possible mass-period correla- 
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tions, we work on the correlation coefficients for the data both in the linear and logarithm 
spaces. A standard method in statistics called the bootstrap method will be used to get the 
confidence intervals of correlation coefficients. 

In the following, we present the construction a new mass-period power-law distribution 
function in Section 2 and the correlated mass-period distribution function in Section 3. 
We present the bootstrap method in Section 4 and describes the results of the correlation 
coefficients and confidence intervals in Section 5. Finally, we provide the concluding remarks 
in the final section. 



2. The Power-Law Distribution Function 

In this section, we construct a new distribution function, assuming that the mass and 
period are two independent power laws. We consider the probability density function (pdf) 
of the power law has the following form 

f P ower(x\k) = C(k)x~ k , 0<a<x<b<oo, (1) 

where the exponent k is an unknown parameter and the constant C(k) is given by the 
normalization requirement that 

rb 



J f power (x\k)dx = ^-(b l ~ k -a l - k ). 



That is, 

Ok) 1 - ' 



b l-k_ a l-k- 

When sampling is from a population described by (1), the parameter k yields the knowl- 
edge of the entire population. Hence, it is natural to seek a method of finding a good es- 
timator of k. The method of maximum likelihood is one of the most popular techniques 
for deriving estimators. Let Xi, ■ ■ • , X n be independent and identically distributed (i.i.d.) 
samples from the pdf f poW er( x \k) , the likelihood function is given by 

n 

L(k> \x\ , • • • , Xn) J^J f power | • 
i=l 

Taking logarithm of both sides, differentiating partially with respect to the parameter k, and 
setting the result to zero we can determine the maximum likelihood estimate (MLE) of k by 
solving 

a 1 - fc -& 1 - fc + (l-Jfe)(& 1 -*ln&-a 1 - fc lna) _ ELi 1 * 1 ^ 



(1 - k)(b 1 ~ k - a 1 -*) n 



(2) 
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Now, we use the power-law distributions 



f^wer( m \km) = C (k m )m km , Cl m < 171 < b m 

and 

fpoweMkp) = C(k p )p- k ", a p <p<b p 

to fit the 233 observed data in the M and P spaces, respectively. First, we choose the range 
of M as follows: 

where M m i n and M max are the smallest and largest mass of the data set. That is, a m = 0.012 
and b m = 18.4. By Eq. (2), we obtain the MLE of k m , which is k m = 0.7805. Similarly, the 
range of P is a p = 1.211909, b p = 4517 .4 and the MLE of k p is k p = 0.9277. 

The histograms of observed data in the M and P spaces are showed in Fig. 1. Fig. 1(a) 
is the one for M and the area covered by this histogram is 115.5. We define the curve, 

fplw^r = 115.5/j^ er (m|fc m ) 

and plot it as the dotted curve in Fig. 1(a) for comparison. Similarly, the histogram in Fig. 
1(b) is for P and its area is 11650. We define another curve, 

fpow P er = ^650fp Ower (p\k p ) 

and plot it as the dotted curve in Fig. 1(b). 



3. The Correlated Distribution Function 

In §2, the distributions of mass and period of the extrasolar planets are described by 
two independent power laws. However, using a data set of 66 exoplanets, Zucker and Mazeh 
(2002) suggested the possible mass-period correlation. Further, our data of 233 samples here 
show that the correlation coefficient in M — P space is 0.1762 and this indicates there exists 
a positive correlation between M and P. It is thus not suitable to use two independent 
power laws to describe the joint mass-period distribution. Therefore, we need to know how 
to simultaneously describe and use probability models to elicit information from the mass 
and period measurements. It is necessary to construct a new distribution function, in which 
the mass and period can be positively correlated and coupled. This was not possible until 
an algorithm for generating two positively correlated beta-distributed random variables was 
provided in Magnussen (2004). This is the reason why we decide to proceed with the beta 
distribution here. 
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The beta distributions are very general, have many possible different functional shapes, 
and have the advantage that the variable boundaries and the normalizations are auto- 
matically considered. The beta distributions are continuous on the finite interval (c, d), 
— oo < c < d < oo, indexed by two positive parameters a, (3. The pdf is given by 

1 (x — c) a-1 (d — x)^ -1 
f be ta(x\a, (3) = ^ _ c y +p ^ , c < x < d, a > 0, (3 > 0, (3) 

where B(a,/3) denotes the beta function, 

B(a,/3) = [ r-^l-tf^dt. 
Jo 

In Eq.(3), the pdf fbeta(x\a, (3) satisfies fbeta(x\a, (3)dx = 1. The beta function B(a, f3) 
can be expressed as 

where the gamma function is 

/•oo 

r(a) = / t a - x e- l dt. (4) 
Jo 

Considering the following transformation 

x — c 

y = i — ' 

a — c 

we then have the pdf 

fXy\a^) = —^—y a -\l- y y-\ 0<y<l, (5) 
B(a,fj) 

which is called the standard beta distribution. 

The beta distribution is often used to model a phenomenon which could be described by 
the values of random variables defined in a finite interval. As the parameters a and (3 vary, 
the beta distribution takes on many shapes, as shown in Fig. 2. The pdf can be strictly 
increasing, strictly decreasing, or U-shaped. The case a — (3 yields a pdf symmetric about 
(c + d)/2. If (3 — 1, the distribution is called a power-function distribution. That is, this 
distribution is one kind of power-law functions. 

From the above discussion, it is clear that the beta distributions are very versatile and 
can be used for many different purposes. This flexibility encourages its empirical use in a wide 
range of applications. For example, Wall et al. (2000) successfully used the beta distribution 
to model both the subgrid-scale pdf and the subgrid-scale Favor pdf of the mixture fraction. 
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Ettoumi et al. (2002) used beta distributions to analyze solar measurements in Algeria. 
Flynn (2004) suggested the beta distribution as a suitable model for human exposure to 
airborne contaminants. Ji et al. (2005) proposed a beta-mixture model to analyze a large 
number of correlation coefficients in bioinformatics. 

The standard beta function is for one variable only. If there are more than one variables 
and they are independent to each other, the extension to multi-variable cases is straightfor- 
ward. However, to have a generalized beta function with a pair of correlated variables is not 
easy. A few algorithms have been proposed to generate pairs of correlated beta-distributed 
random variables numerically (please see Johnson 1987, Loukas 1984, Michael & Schucany 
2002). Due to the limitations of these algorithms, they can only be used for particular types 
of data set. The algorithm in Magnussen (2004) can generate a pair of positively correlated 
beta-distributed random variables without any limitations. We thus use it to construct the 
numerical mass-period distribution function f(M,P\a m ,j3 m ,a p ,/3 p ). 

The probability that a planet with mass and orbital period in the range [M, M + dM], 
[P, P + dP] is given by 

f(M, P\a m , (3 m , a p , (3 p )dMdP, 

where the marginal distributions of M and P follow the beta distribution as in Eq. (3) with 
parameters (a m ,/3 m ) and (a p ,fl p ), respectively. From the data, the boundaries can be set 
such that mi < M < m 2 , p\ < P < p 2 . Then, the marginal distributions of the random 
variables 

Mi = — - — and PI = £Z£1 
m 2 - mi p 2 - pi 

should follow the standard beta distributions as Eq.(5) with parameters (a m , j3 m ) and (a p , (3 P ), 

respectively. 



Let us define that 



and 



where 



<Si = p(Mi, Pi) x (1 + a m + a p ) x C (6) 
5 2 = p(M 1 , Px) x (1 + f3 m + p p ) x C, (7) 



c V a m ap(3m(3 p (l + a m + (3 m )(l + a p + (3 P ) 

(1 + a p )(l + (3 m )(l + (3 P ) + a m (l + (3 m + (3 p + (3 m (3 p + a p (l + (3 m + (3 P )) ' 

and p(Mi,Pi) is the correlation coefficient between Mi and Pi. Then, M 1 and Pi generated 
by below equations would be a pair of correlated beta-distributed variables: 

Mi = OKJ + GW (9) 
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and 

G(a* p ) + 

Pl = G(a;) + G(5 1 ) + G(P;) + G(5 2 y (10) 

where G(a) is a random variable distributed as a gamma distribution with parameters a 
and 1 and a* m) (3^, a*, (3* are defined by 

(11) 
(12) 
(13) 
(14) 

Note that the pdf of a gamma distribution with parameters a and (3 is (Hogg & Craig 1989) 

f( x ) = — i- x^e-'W, (15) 

where < a; < oo. 

The above procedure to generate a pair of positively correlated beta-distributed M x and 
P 1 variables can be summarized as: 

Magnussen Algorithm 



„ * 




-5i 


13* 


= !3 m 


-5 2 


„ * 


= ot p ■ 


-6! 


13* 


= P P - 


-5 2 . 



Step 1 Assume that the marginal distribution of Mi, i.e. fM 1 (m\a m , p m ), is a standard beta 
distribution with parameters (cx m ,(3 m ). Through the maximum likelihood method, we 
employ the data to get the best estimation (a m , /3 m ) of (a m ,/3 m ). Similarly, we also 
get the best estimation (a p , /3 P ) of (a p , (3 P ) for P 1 . 

Step 2 Calculate the value of C by Eq.(8). 

Step 3 Calculate the correlation coefficient p(M±, Pi) from the data and use it as the value of 
p(M 1 ,P 1 ). 

Step 4 Calculate S 1: S 2 by Eqs.(6)-(7). 

Step 5 Calculate <, (3* m , a* p , (3* p by Eqs.(ll)-(14). 

Step 6 Generate pairs of M x , Pi by Eqs.(9)-(10). 
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We now apply the above algorithm on the data set of 233 exoplanets. To avoid 
the possible singularity, the range of M is chosen to be (M min /1.5,1.5M max ). That is, 
m\ = 0.008, m 2 = 26.7. By the same reason, the range of P is taken as pi = P m j n /1.5 = 
0.8079, £>2 = l-5Pmax = 6776.1, where P m in and P max are the smallest and largest of period 
of the observed data. Thus, the MLE of a m , (3 m , a p and (3 P are a m = 0.6524, f3 m = 5.9076, 
dip = 0.3697, P P = 3.8445, respectively. In addition, the data shows that mass-period corre- 
lation coefficient p(Mi,Pi) = 0.1762. We then get all the parameters' values as C = 0.2092, 
5i = 0.0745, 5 2 = 0.3961, and a* m = 0.5779,/?^ = 5.5115, a* = 0.2952, (3* p = 3.4484. 

Because the area of the histogram for M in Fig. 1(a) is 115.5, we define the curve, 

fZZ? = 115.5 f beta (m\a m J m ) 

and plot it as the solid curve in Fig. 1(a) for comparison. Similarly, the area of the histogram 
for P in Fig. 1(b) is 11650, so we define another curve, 

C: = 11650 f beta (p\a p Jp) 

and plot it as the solid curve in Fig. 1(b). 

These plots indicate that the beta distribution presents a better fitting with the data, 
comparing with the power law. Due to that there is no closed form for the positively corre- 
lated beta distribution, we numerically plot the joint distribution of M and P as shown in 
Fig. 3. Fig. 3 presents the three dimensional plot of our correlated mass-period distribution 
function. Fig. 4 shows the contour of it in smaller ranges of mass and period. We thus 
successfully construct, for the first time in this field, the correlated mass-period distribution 
function. 

Please note that, for pairs of quantities (xi,yi),i = 1, . . . , n, the correlation coefficient 
9 is usually given by 

9= - ^t^-^-y) (16) 

where x = Yli=i x il n i V = Y^=i yi/ n - The value of 9 lies between —1 and 1. If the data points 
completely lie on a straight line with positive (negative) slope, the correlation coefficient 9 = 
1 {6 = — 1) and it is called "complete positive correlation" ("complete negative correlation"). 
When the data points are randomly distributed, the variables x and y are uncorrected and 9 
is near zero. Thus, the value 9 is regarded as one conventional way to quantitatively describe 
the strength of relationship between x and y. Thus, our p(Mi,Pi) = 0.1762 was obtained 
by the above equation. 
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4. The Bootstrap Confidence Intervals 

To assess the statistical significance of the possible correlation, it would be good to 
have the confidence interval corresponding to a given confidence level. We use the bootstrap 
method to construct confidence intervals here. Statistically, we determine the characteristics 
of the population by taking samples. Since the sample represents the population, analogous 
characteristics of the sample should give us information about the population characteristics. 
The bootstrap method proposed by Efron (1979) is a simple and straight-forward method to 
calculate the approximated biases, standard deviations, and confidence intervals, for exam- 
ple. It gives the population characteristics by taking samples repeatedly from the original 
data set. 

The bootstrap method is used for i.i.d. data. DiCiccio & Efron (1996) found that the 
bootstrap confidence intervals are more accurate than the classical normal approximation 
intervals. The standard bootstrap method for confidence interval can be described as follows. 
Given an observed i.i.d. sample z = {zi, • • • , z n } from an unknown distribution function F, 
we want to construct a confidence interval for an interesting parameter 9 = 9(F) based on z. 
Let F be the empirical distribution function which is defined to be the discrete distribution 
that assigns the probability 1/n on each value Zi, i — 1, • • • , n. The key idea of the bootstrap 
method is a bootstrap sample, which is defined to be a random sample of size n drawn from 
F, say z* = {z{ , • • • , z^}. That is, the bootstrap data points z*, ■ ■ ■ , z* n are a random sample 
of size n drawn from the data set {z±, • • • , z n }. 

Let 9 = s(z) be an estimate of 9. Corresponding to a bootstrap data set z* is a bootstrap 
replication of 9, 

9* = s(z*). 

The quantity s(z*) is the result of applying the same function s(-) to z* as was applied to 
z. For example, if s(z) is the sample mean z = *YTi=\ z il n then s(z*) is the mean of the 
bootstrap data set, z* = Y^i=i z i l n - The bootstrap algorithm, described next, is a data- 
based simulation procedure to obtain a good approximation of the confidence interval for 
9. 

Bootstrap Algorithm 

Step 1 Draw a " bootstrap sample" z* = (z*, • • • , z*) according to F. 

Step 2 Evaluate 9* = s(z*), where s(z*) is the value of s(z) based on z*. 

Repeat the previous two steps a large number of times, say B times, to obtain 0*, • • • , 9* B . 

Step 3 Sort 9\, • ■ • , 9* B to be the ordered list 9* {l) , ■ ■ ■ , 9f B) . 
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Step 4 Let 9* B a \ < a < 0.5, be the (100a)th empirical percentile of the b 
that is, the (a_B)th value in the ordered list of 6%\, • • • , 9(b)- 
Likewise, let 9* B a ^ be the the (1 — a)Bth value in the ordered list of 9*^ 

Then, the approximate 1 — 2a confidence interval for 9 is (9* B a \ 9* B l a ^). 

Moreover, if Ba is not an integer, the following procedure can be used. Assuming 

< a < 0.5, let k = [(B + l)a], i.e. the largest integer < (B + l)a. Then we define the 9* B 
and 9* B {1 ~ a) by the fcth and (B + l- k)th values of 9* {b) , b = 1, • • • , B. 

In this paper, the bootstrap algorithm is used to construct confidence intervals for 
correlation coefficients. Assuming that n objects are sampled from a population and two 
numerical characteristics are measured on each of them, we end up with bivariate random 
sample points z\ = (x%, yi), ■ • ■ , z n — (x n , y n ). Let 9 be the population correlation coefficient. 
Based on n data points z%, ■ ■ ■ , z n , the population correlation coefficient 9 is estimated by the 
sample correlation coefficient 9 defined in Eq.( ll6jl . Independent repetitions of the bootstrap 
sampling process give B bootstrap replications 8*, ■ ■ ■ , 9* B . Then we obtain the approximate 

1 — 2a confidence interval for 9, which is (9* B a \ 9* B l a ^). According to Efron & Tibshirani 
(1993), the number of independent repetitions of the bootstrap sampling shall be set as 
B = 2000. 



5. Mass-Period Correlations 

We first work on the mass-period correlation inM-P space. With the 233 observed 
data, we calculate the correlation coefficients and use the bootstrap method to determine 
the confidence intervals. We find that the correlation coefficient between M and P is 0.1762 
and the resulting 95% bootstrap confidence interval is (0.0575,0.3130). This indicates that 
the mass M and period P has a weak positive correlation. Fig. 5 shows the 233 samples in 
M — P space. In deed, it looks as that the distribution is not completely random. However, 
the positive correlation is difficult to recognize and so it is consistent with a weak correlation. 

On the other hand, for the mass and period distributions in InM-lnP space, the cor- 
relation coefficient between InM and InP is 0.3876 and the corresponding 95% bootstrap 
confidence interval is (0.2668, 0.5001). It clearly indicates that there is a positive mass-period 
correlation for the data in In M-ln P space. 

Fig. 6 shows the 233 exoplanets in In M-ln P space and it seems there are many fine 
structures of data distributions. For example, for those points with InM > 0, it is far more 
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crowded in the region with InP > 5. This is, of course, related to the deficit of massive 
close-in planets, which makes the main contribution on the positive mass-period correlation. 
On the other hand, partially due to the effort of transit surveys, the discovery of many hot 
Jupiters also makes it a bit crowded around (InM, InP) = (—0.5, 1). However, the deficit of 
sub- Jupiter mass planets at separations about 0.5 AU mentioned in Papaloizou & Terquem 
(2006) seems to disappear in this plot of 233 data points. 

In fact, the mass-period correlation was theoretically studied in Patzold & Rauer (2002) 
and Jiang et al.(2003), which focused on the paucity of massive close-in planets. The ex- 
planation for this paucity and the correlation could be that the tidal interactions with host 
stars make the massive close-in planets migrate inward and finally merge with the stars. 

6. Concluding Remarks 

In this paper, we first construct a new mass-period distribution function of exoplanets, 
in which the correlation is considered. This is the first time in this field and was not possible 
until the method was proposed in Magnussen (2004). 

The correlation coefficients of exoplanet data in M-P space and In M-ln P space are 
further determined and the bootstrap method is then used to construct the confidence in- 
tervals of correlation coefficients at particular confidence levels. We confirm that there is a 
mass-period correlation for exoplanets. In addition to the paucity of massive close-in planets, 
there are other fine structures in the data distribution to be investigated in the future. 
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Fig. 1. — The mass and period histograms: (a) the mass histogram of 233 exoplanets, where 
the solid curve is an d the dotted curve is fp^ r ; (b) the period histogram of 233 

exoplanets, where the solid curve is flll^ and the dotted curve is fp^f er . Please note that 
the mass's unit is Jupiter Mass (Mj) and the period's unit is days. The bin size of the above 
mass histogram is 0.5 Mj and the bin size of the above period histogram is 50 days. 
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Fig. 2. — Beta distribution functions in Eq.(3) with c = 1, d — 5: the solid curve is 
for a = (3 = 0.5; the dotted curve is for a = 0.5,/? = 2; the short dashed curve is for 
a — 2, (5 — 0.5; and the long dashed curve is for a = f3 = 2. 
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Fig. 3. — The correlated mass-period distribution function. The unit of the mass M is 
Jupiter Mass (Mj) and the unit of the period P is days. 
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Fig. 4. — The contour of the correlated mass-period distribution function. The unit of the 
mass M is Jupiter Mass (Mj) and the unit of the period P is days. 
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Fig. 5. — The data of 233 exoplanets in M — P space. The unit of the mass M is Jupiter 
Mass (Mj) and the unit of the period P is days. 
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Fig. 6. — The data of 233 exoplanets in InM — InP space. The unit of the mass M is Jupiter 
Mass (Mj) and the unit of the period P is days. 



